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Abstract. In this paper we consider Galerkin-finite element methods that approximate the solutions of 
initial-boundary-value problems in one space dimension for parabolic and Schrodinger evolution equations 
with dynamical boundary conditions. Error estimates of optimal rates of convergence in and are 
proved for the accociated semidiscrete and fully discrete Crank-Nicolson-Galerkin approximations. The 
problem involving the Schrodinger equation is motivated by considering the standard 'parabolic' (paraxial) 
approximation to the Helmholtz equation, used in underwater acoustics to model long-range sound propa- 
gation in the sea, in the specific case of a domain with a rigid bottom of variable topography. This model 
is contrasted with alternative ones that avoid the dynamical bottom boundary condition and are shown 
to yield qualitatively better approximations. In the (real) parabolic case, numerical approximations are 
considered for dynamical boundary conditions of reactive and dissipative type. 

1. Introduction 

Our main goal in this paper is to analyze Galerkin-finite element methods for initial-boundary-value 
problems, involving dynamical boundary conditions, for the linear Schrodinger and the heat equations. In 
addition, in a specific problem arising in underwater acoustics and modelled by the Schrodinger equation, 
we will also consider an alternative boundary condition and evaluate, analytically and numerically, the two 
models. 

We start with the underwater acoustic application. Consider the Helmholtz equation (HE) in cylindrical 
coordinates in the presence of cylindrical symmetry 

(HE) Ap + A:gr;2(r,z)p = 0. 

Here z > is the depth variable increasing downwards and r > is the horizontal distance (range) from a 
harmonic point source of frequency /o placed on the z axis. For simplicity we shall assume that the medium 
consists of a single layer of water of constant density, occupying the region, < z < ^(r), r > 0, between 
the free surface z = and the range-dependent bottom z = £{r) (see Fig. 1); £ — £{r) will be assumed to 
be smooth and positive. The function p = p{r, z) is the acoustic pressure, ko — is a reference wave 
number, co a reference sound speed, and ri{r, z) the index of refraction, defined as j^^i where c(r, z) is the 
speed of sound in the water. (HE) is supplemented by the surface 'pressure-release' condition p(r, 0) = 0. In 
the case of a soft bottom the homogeneous Dirichlet boundary condition 

(D) p = at z = £{r) 

is assumed to hold. The case of a rigid bottom is modelled by a Neumann boundary condition (with £ = 
(N) p^ - £{r)pr =0 at z = £{r). 
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Applying the change of variables p{r, z) — ipir, z) , assuming that |2iA:o'0r |^|'0rr| (narrow-angle paraxial 
approximation) and neglecting terms of O(^) (far-field approximation) we arrive (cf., e.g., [22]; [H], [S]) at 
the standard 'Parabolic' Equation (PE), which is a linear Schrodinger equation of the form 

(PE) Vr= 2k 

where ip = ilj{r, z) is a complex- valued function of the two real variables r and z. The (PE) has been widely 
used in underwater acoustics to model one-way, long-range sound propagation near the horizontal plane of 
the source, in inhomogeneous, weakly range-dependent marine environments. Its solution will be sought in 
the domain < z < £{r), r > 0. The (PE) will be supplemented by an initial condition tp{Q,z) = tpo{z), 
< z < £(0), modelling the source at r = 0, the surface boundary condition = foi z — 0, r > 0, and 
a bottom boundary condition obtained by transforming (D) or (N). The Dirichlet boundary condition (D) 
remains of the same type (V' = at z = ((f)) while the Neumann boimdary condition (N) is transformed to 
a condition of the form 

(PN) Tjj,-i{r)^Pr-9Bir)eir)^P^0 at z^e{r), 

where gsif) is complex- valued and is usually taken simply as i/co- 

The theory and numerical analysis of this initial-boundary- value problem (ibvp) with the Dirichlet bottom 
boundary condition is standard, cf.e.g. [19] . [3], and will not be considered any further. On the other hand 
the analysis is complicated in the case of the Neumann boundary condition, when £(r) is not the zero 
function, due to the presence of the term ipr in (PN). In [1] Abrahamsson and Kreiss proved existence and 
uniqueness of solutions, in the case of a strictly monotone bottom, i.e. when i(r) is of one sign for r > 0. 
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Figure 1. The domain of the initial boundary value problems for the (PE) in the r, z variables 



We shall transform the above ibvp's to equivalent ones posed on a horizontal strip. With this aim in mind, 
we first introduce non-dimensional variables as in [4], defined hy y := j^, t := w := where we take 

L := ^ and Vrcf max|i^o|- Then, letting s{t) fco£(^), git) = fco5i,(^), l{t,y) fl^'l^'i) ^ l] 
we see that the (PE) becomes 

(1.1) wt = ^Wyy + ij{t,y)w, 0<y<s{t), t>0. 

We note that the index of refraction ?/, and consequently the function 7, may be taken to be complex-valued 
in order to model attenuation of sound in the water. The initial condition becomes 

(1-2) w{0,y)=wo{y):^-^Mi;)^ < y < s(0). 

The surface condition remains the same, i.e., 

(1.3) w{t,0)=0, <>0, 
while the boundary condition (PN) becomes 

(1.4) wy{t, s{t)) - s{t) [wt{t, s{t)) + g{t) wit, s{t))\ =0, t>0. 

We now perform the range-dependent change of depth variable x := -j^, that maps the domain of the 
problem onto the horizontal strip 0<a;<l,<>0. We also make the transformation 

(1.5) u{t,x) = exp{—iS{t)x'^) w{t, s{t)x), 



which defines the new field variable u{t,x) for < a; < 1, t > 0. In dLS]) S{t) := ii^Ml^ t > 0, where a dot 
denotes differentiation with respect to t. In terms of the new variables (|1.1[) becomes 

(1.6) ut = ia{t)u^^ +il3{t,x)u, < x < 1, t > 0, 

where, for < x < 1, i > 0, 

'*(^) = 25%)' l3{t,x) = f3n{t,x) +iPj{t,x), 
^^■^^ Mt,x)^Re[^{t,xs{t))]~^'^ii^, I3,{t,x)=lm[j{t,xsm + i^- 

The purpose of introducing in (|1.5p the factor with 5 = ^ is to avoid the presence of a Ux term in 

the right-hand side of (|1.6[) and, consequently, simplify somewhat the analysis. Under the transformation 
()1.5|) . the initial and boundary conditions (|1.2p - (|1.4p change accordingly. Specifically, we have 



(1.11) s,{t) :^ ^ \\)' , m{t):^git)s,it)+iis,{t)Sit)-2S{t)), < > 0. 



(1.8) u(0,a;) =Mo(a;) :=e-''^(°)^'wo(2;s(0)) Vxe[0,l], 

(1.9) M(i,0) = 0, i>0, 
and 

(1.10) Uxit,l) ^ si{t)ut{t,l) +m{t)u{t,l), i>0, 
where 

1 + S(t)2 

The boundary condition (jl.lOp is an example of a dynamical boundary condition, because it involves (if 
s ^ 0) the value of ut at the boundary. As was already mentioned, the well-posedness of ibvp's of the type 
{ (|1.6p . (jl.Sp . (|1.9p . (|1.10p }. for t in a finite interval [0,T], was proved in [IJ under the assumption that s{t) 
is of one sign for all t G [0,T]. One of our main purposes in this paper is to construct and analyze fully 
discrete Galerkin-finitc element methods for approximating the solution of the above ibvp. 

We consider the ibvp consisting of (|1.6p - (|l.lip . We assume that the bottom is upsloping, i.e. that s{t) < 0, 
and that the problem has a unique solution, smooth enough for the purposes of the error estimation. In 
paragraphs 2.1 and 2.2 we discretize the problem in x by the standard Galerkin method and prove optimal- 
order and estimates for the error of the resulting semidiscretization. This is achieved by using 
appropriate properties of the and the elliptic projections onto the finite element subspace and a relevant 
superconvergence result. (The difficulty of the problem lies in the presence of the ut term in (|1.10p : the 
condition s{t) < 0, which implies that si{t) < 0, is needed to obtain a basic energy inequality for the error 
of the semidiscretization). Subsequently, in paragraph 2.2, we discretize the semidiscrete problem in the t 
variable using a Crank-Nicolson type method with a variable step-length. Again, under the assumption that 
s{t) < for < t < T, we prove and error estimates which are of optimal order in x and t. 

In order to overcome the analytical and numerical difficulties caused by dynamical boundary conditions 
of the form p.lOp Abrahamsson and Kreiss proposed in |l2j an alternative rigid bottom boundary condition, 
which, in the case of (PE), is of the form 

(AK) V^-ifcoi(r)V' = at z ^ £{r). 

This condition may be viewed as a 'paraxialization' of (PN). When the nondimensionalization z y, r —> t, 
"0 ^ w is performed, (AK) becomes 

(1.12) Wy{t,s{t)) -is{t)w{t,s{t)) ^0, t>0. 

Finally, after changing the depth variable by x — y/s{t) and the dependent variable by (jl.Sp . it is not hard 
to see that (11.121) becomes simply 

(1.13) Ux{t,l)^Q, t>0. 

The proof of the well-posedness of the ibvp consisting of (|1.6p - (|1.9l) and p.l3p is standard, cf. Its 
numerical analysis too is straightforward; under no restriction on the sign of s{t) we prove in paragraph 
2.3 optimal-order and error estimates for the standard semidiscrete Galerkin scheme and its Crank- 
Nicolson full discretization. 
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In Section[3]we present results of various numerical experiments that we performed for problems on variable 
domains with the Neumann and Abrahamsson-Kreiss bottom boundary conditions, using the fully discrete 
finite element methods analyzed in Section [21 As predicted by the theoretical stability and convergence 
analysis, the finite element scheme is stable and second-order accurate when Neumann boundary conditions 
are considered in domains with upsloping bottoms. (It also appears to be convergent in small scale problems 
with downsloping bottoms and also in more realistic examples if the downsloping bottom has very small 
slope.) The scheme with the Abrahamsson-Kreiss condition behaved well, as predicted by the theory, in all 
examples of bottoms of arbitrary shape that we ran. 

When we compared the results of the schemes using both boundary conditions in the case of the upsloping 
and downsloping rigid bottom ASA wedge (a standard test problem for long range sound propagation in 
underwater acoustics, [15j). we found that in the upsloping case there was very good agreement between the 
two schemes. In the downsloping case, the scheme implementing the Neumann boundary condition was not 
convergent. This is in agreement with the results of Abrahamsson and Kreiss, [T], [2], who pointed out that 
for some downsloping bottom profiles one may observe instabilities in the case of the Neumann boundary 
condition. On the other hand, the scheme with the Abrahamsson-Kreiss condition was convergent and its 
results agreed well with those furnished by the finite difference code IFD, ^l6], [17], [H], implemented with 
the rigid bottom boundary condition option. The IFD scheme uses a discretized version of the Neumann 
boundary condition (PN), wherein the 4'r term is replaced by the right-hand side of the (PE). We prove a 
priori estimates for the resulting ibvp. 

A final point of interest emerging from the numerical experiments is that, for some downsloping bottom 
profiles s{t) with an inflection point at some t — t* , we observed violent growth of the L^-norm of the 
numerical solution of the problem with the Neumann boundary condition for t > t* . This growth (blow-up?) 
of the solution seems to be a feature of the problem and not an artifact of the numerical scheme. 

Error estimates for a finite difference scheme of second-order of accuracy in x and t for some of the 
ibvp's considered here were proved in In the case of the Neumann boundary condition (jl.lOp these error 
estimates were shown to hold not only when s{t) < but also in the strictly downsloping case s{t) > 0, 
t G [0, T], as a result of the validity of a certain discrete estimate; this estimate mimics an analogous 
estimate for the continuous problem, which holds provided s{t) < or s{t) > when t g [0, T]. In [21] Sturm 
considered the Abrahamsson-Kreiss condition for the (PE) in three dimensions over a variable bottom in the 
more general case of a multilayered fluid medium with homothetic layers. When restricted to single layer 
problems in the presence of azimuthal symmetry, the scheme of [21 is similar to the one analyzed here in the 
case of the Abrahamsson-Kreiss bottom boundary condition. We have considerably modified the analysis 
of [21' and obtain optimal-order estimates, since, by using the transformation ()1.5|) . we essentially avoid an 
elliptic projection with time-dependent terms. We finally mention that a uniform range step version of the 
scheme of this paper and also three-dimensional extensions thereof were analyzed in [5] . 

The problem addressed in the present paper, namely sound propagation modelled by the (PE) in a single 
layer of water over a rigid bottom, is, of course, an idealized model problem in underwater acoustics. More 
realistic environments consist, for example, of a layer of water above several layers of fluid sediments of 
different density, speed of sound and attenuation overlying a rigid or soft bottom. If the layers are separated 
by interfaces of weakly range-dependent topography and low backscatter is expected, long-range sound 
propagation may again be modelled by the (PE) in each layer with transmission conditions (continuity of 
ip and of where p is the density and n the normal direction to the interface) imposed across the layer 
interfaces. Hence, the issue arises of how to treat the dynamical interface condition, now involving ipr on 
both sides of an interface, and the ensuing problems are analogous to those encountered in the case of the 
dynamical bottom boundary condition. The analysis is more complicated now, as it appears that possible 
non-homotheticity of the layers has to be balanced by the jump across the interface in the imaginary part of 
the analog of the function 7, cf. (|l.ip . in order to ensure the well-posedeness of the problem, [T3]. For a recent 
review of several issues regarding the interface problem for the (PE), we refer the reader to [12j : references 
to underwater acoustics computations with the (PE) in the presence of interfaces with change-of-variable 
techniques include e.g. |4], [21] and [12]. Here we just wish to point out that range-dependent topography has 
often been approximated in practice by 'staircase' (piecewise horizontal) bottoms and interfaces. This raises 
the issue of what boundary / interface conditions to pose on the vertical part of the steps of the staircase. 
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Moreover, it is well documented that staircase approximations lead to nonphysical energy losses or gains, cf. 
e.g. [TS], [50]. To alleviate this problem of energy non-conservation, change-of- variable techniques may be 
used as in the present paper. They may also be extended to interface, JJj, or 3D-problems, pT| . 

We turn now to one-dimensional (real) parabolic problems with dynamical boundary conditions. We 
consider the following model problem: For < T < cx) we seek a real-valued function u — u{t, x) defined for 
{t, x) e [0, T] X [0, 1] and satisfying 



where a{t) > a* > for t e [0,T] and /, e, 6, g, uq are smooth, real-valued functions. Such problems 
occur in heat conduction, ([llj. Section 4.3.5), and in other areas; see [14] for a fuller list of references. Our 
aim is to construct fully discrete Galerkin-finite element approximations for the ibvp (|1.14p and prove error 
estimates, with techniques analogous to those used in the case of the Schrodinger equation. We consider two 
different cases depending on the sign of the function e in the dynamical boundary condition. 

We treat first the dissipative case, characterized by the hypothesis that e{t) < for all t S [0,T], and in 
which the ibvp (|1.14[) is well posed, cf. e.g. [14]. In paragraph 4.1, applying the standard Galerkin method 
to this case, we prove optimal-order and estimates for the error of the resulting semidiscretization 
and for the Crank-Nicolson-Galerkin fully discrete scheme. Matters are more complicated in the reactive 
case, wherein e{t) > for t G [0, T]. In this case the problem is well posed in one space dimension as in the 
case at hand, but in general is not well posed in higher dimensions, [21], [7|. To construct a Galerkin-finite 
element method in this case, we replace the term ut in the dynamical boundary condition using the p.d.e. 
in (|1.14p . thus obtaining a boundary condition involving Uxxit, 1)- The resulting ibvp is discretized in space 
by means of a H^-type Galerkin method that uses finite element spaces consisting of piecewise polynomial 
functions in of degree at least three. In paragraph 4.2 we analyze this method and prove optimal-order 
error estimates for the semidiscrete approximation and the fully discrete one when the Crank-Nicolson 
scheme is used in time-stepping. The case where e(t) changes sign in [0, T] is under investigation; for a 
discussion cf. [8]. 



2.1. Preliminaries. Let D := (0,1). We will denote by L^(D) the space of the Lebesgue measurable 
complex-valued functions which are square integrable on D, and by || • || the standard norm of L'^{D), i.e., 
11/11 := {/^ \f{x)\'^dx}^ for / e L'^{D). The inner product in L'^{D) that induces the norm || • || will be 
denoted by (•,•), i.e. (/i,/2) := fi{x)J^ dx for /i, /a G L^{D). Also, we will denote by L°°{D) the 
space of the Lebesgue measurable functions which are bounded a.e. on D, and by | • |oo the associated 
norm, i.e., |/|oo := esssup^^ |/| for / e L°°{D). For s G No, we denote by H'^{D) the Sobolev space of 
complex- valued functions having generalized derivatives up to order s in L'^{D), and by || • ||s its usual norm, 

i-e- 11/11. := {ELo W^iffY for / G H'{D). In addition, we set \v\i := for v e H^{D). Also, Hi(L>) 
will denote the subspace of H^{D) consisting of functions which vanish at cc = in the sense of trace; we 
set WiD) = H^iD) n W^^D) for s > 2. In addition, for s G No, we denote by VF^'°°(D) the Sobolev space 
of complex-valued functions having generalized derivatives up to order s in L°°{D), and by | • |s,oo its usual 
norm, i.e. |/|s,oo := inaxo<^<s |9^/|oo for / G W^'°°{D). In what follows, C will denote a generic constant 
independent of the discretization parameters and having in general different values at any two different 
places. 

For later use, we recall the well-known Poincare-Friedrichs inequality 



(1.14) 



Ut = a{t)uxx + f3{t, x)u + f{t, x) V(t, x) G [0, T] x [0, 1] 
u{t,0) = yte [0,T], 

a(t)u,At, 1) = e{t)ut{t, 1) + 6{t)u{t, 1) -I- g{t) Vi G [0, T 
u{0,x) — uo{x) Vx G [0,1], 



2. Numerical Schemes and Error Estimates for the (PE) 




ll^ll < Cp^ |f|i G m\D), 



\v\o. < \v\i "^veM^D) 
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and the trace inequality 

(2.3) \v{l)\^ < 2\\v\\ \v\i Vv e m\D). 

Let r e N and Sh be a finite dimensional subspace of (D) consisting of complex- valued functions that 
are polynomials of degree less or equal to r in each interval of a non-uniform partition of D with maximum 
length h € (0, /i^,]. It is well-known, [10] . that the following approximation property holds: 

(2.4) M {\\v~x\\+h\\v-x\\i} <Ch'+^\\v\\,+u y v € U'+\D), y h e {0, K], s = 0, . . . , r. 

Also, we assume that the following inverse inequality holds 

(2.5) < C h-^ \m Vcf^Sh, yhe{0,h,], 

which is true when, for example, the partition of D is quasi-uniform, [lOj . In addition, we define the 
i^— projection operator : L'^{D) Sh. by 

iPhV,^) = iv,cl)) yc^eSh, yveL^iD), 

and the elliptic projection operator Rh : H^{D) Sh by 

(2.6) B{RhV,c^)^B{v,c^) y<PeSh, yveH^D), 

where B is the sesquilinear form defined for u, w E H^{D) by B{u, w) := {u' , w'). It follows, [TU], [23]) that 

(2.7) \\RhV-v\\+h\\RhV-v\\i< Ch'+^\\v\\s+i, yvem'+\D),yhe{0,h^], s = 0, . . . , r. 
Finally, for v E L^{D), we define the discrete negative norm 

:-sup{M : ^QSh and 0^o}, V/iG(0,/i,]. 
Lemma 2.1. The elliptic projection operator Rh has the following property: 

(2.8) Rhv{i)^v{i), yvem^{D). 

Proof. Let v e M^{D) and lu be the element of 5*^, given by u!{x) = x for x E D. Then (|2.6p gives 
Rhv{l) — v{l) = B{RhV — v,uj) = 0, which is the desired result. □ 

Lemma 2.2. Let uj e C^D). Then 

(2.9) |Ph(^0)|i < C|w|i,oo|0|i y^beSh, V/iG(0,M- 



Proof Let h G (0, /i*] and (j) G Sh- Since |Pft(w0)|i < \Ph{L0<j)- Rh{uJ(t)))\i + |i?/»(w(^)|i, using ([211) and (12^]) 
we arrive at \Ph{u!(t))\i < C h^^ \\^(l>— -R/i(w0)|| -|- Iwc/)!!. Next, we use the estimate (|2.7p for s = to obtain 
\Ph{u}4>)\i < C [\uj\oo \(l>\i + |w'|oo||0|| ]. Thus, the bound (|2.9p follows by combining the latter inequality 
and (Hi]). □ 

2.2. The Neumann (dynamical) boundary condition. In this subsection, we shall consider the (PE) 
with the Neumann boundary condition, i.e. the ibvp (|1.6|) . (II. 8p . (|1.9p . p.lOp . We shall write this problem 
in a slightly more general form, as follows. For T > given, we seek a function u : [0,T] x D ^ C satisfying 

ut = i a(t) + i I3{t, x)u + fit, x) V (i, x) G [0, T] x D, 

u{t,0) = yte [0,T], 

u^{t,l)^fi{t)[S{t)ut{t,l) + G{t)u{t,l)] yte[Q,T], 
u{0,x) = uo{x) Vx G D. 

We shall assume that a : [0,T] M\{0}, (3, f : [0,T] x D ^ C, uq : D ^ C, fi, S : [0,T] ^ R and G : 
[0, T] ^ C are given functions. We shall assume that the solution u of (Af) exists uniquely, and that the data 
and the solution of (A^) are smooth enough for the purposes of the error estimates that will follow. (In some 
numerical experiments of Section [3] we shall revert to the specific physical data in (|1.9p , (jl.lOp , (jl.lip , and 
take the functions a(t), /3(i,x) as in ([LHI), ^J.{t) = S{t) = ttW)7' '^^^^ ^ 9(t) S{t) + i[S{t) S{t) - s'^{t)], 
where (5 = ^•) 



2.2.1. Semidiscrete approximation. The weak formulation of (A/"), obtained by taking the L^{D) inner prod- 
uct of the p.d.e. in (A/") with a function in H^(Z?), integrating by parts and using the boundary conditions, 
motivates defining Uh '■ [0,T] — > Sh, the semidiscrete approximation of m, by the equation 



{dtUhit, •), 4>) = iaW A^W [ S{t) dtu^it, 1) + G(t) u^it, 1) ] 0(1) 

-ia(i)6(u^(i,-),0)+i(/3(^,-)^/.(i,-),0) + (/(i, •),'/') V0e^,„ Vie[0,T], 

and 

(2.11) u,,(0,-) = i?hwo(-)- 

Proposition 2.3. The problem (|2.10p - (|2.1ip admits a unique solution Uh E C^{[0,T]; Sh)- 

Proof. Let dim(S'/i) — J and {'/>j}j^i be a basis of 5/i consisting of real-vahied functions. Hence, we 
have RhUo = Y.]=il°j 4>j and Uh{t,x) = Y.''j=ilji't) 4>j{x), where 7j : [0,T] ^ C for j = 1,...,J. 
Then ([2T0| - (l2ll|) is equivalent Jo the foUowmg o.d^e. initial-value problem:^ Find G G Ci([0,r];CO 
such that G(0) = G"^and I(t)G'(t) = B{t)G{t) + F{t), yt e [0,r], where G{t) := (7i(t), . . . , 7^,(t))^, 
GO :=j70,...,7?)^, A : [0,T] ^ C'-' with A,,(t) := (</),,</.,)- ia(t) 5(t) A*(i) '/'^l) '/'^(l), B : [0, T] C^^ 
with Bej{t) := ia{t) G{t) ^i{t) (t)e{l) (jj^il) - ia{t) B{<j)i, (jj-j) + i{P{t, ■) <j)i, (t)j), and F : [0,T] ^ C'" with 
F{t) := ((/(i, •), (j>i), . . ■ , {f{t, •), cfij))'^ ■ Since A, B, are continuous maps, to ensure existence and unique- 
ness of the solution G, it is sufficient to show that A{t) is nonsingular for t G [0, T]. Indeed, letting t e [0, T] 
and X S Kei{A{t)), we have Ke{x'^ A{t)x) = 0, from which we conclude that || X]j=i ^j'^iH — ^ and hence 
a; = 0. □ 

Let us first present a superconvergence error estimate for the semidiscrete approximation Uh. 

Proposition 2.4. Let u be the solution of (J\f) and Uh its semidiscrete approximation defined by (j2.10p - 
(f2TT|) . Assume that ^{t) < and S{t) > for t e [0, T]. Then 

(2.12) \\uhit,-)^RHu{t,-)\\i< Gh^+^ (^J^T^{T)dTy Vte[0,T], V/ie(0,/i.], 
where T^{t) := ||u(r, •)||2^i + \\dtu{T, ^H^+i + ELo Io W^t^is, Ollr+i o^s- 

Proo/. Let h e (0, /i*]. Oh := it/i - and ^(t) := Using and ([211) we obtain 



{dtehit,-),^)^iait)^iit) [S{t)dteh{t,l) + Git) 9hit,l)] 0(1) 

(2.13) - ia(t) 6(0^ (t, •), 0) + i [Phm, •) ^^/.(^, •)), ^) 

+ (*4t,-),0) V0e5ft, ViG[0,T], 

where ^I'* := — i?;i(9fu)] —if3{u — Rhu). Set = dt9h in (|2.13p and then take imaginary parts to obtain 

jr,\On{t, < \^l{t)\ [-25* \dMt, 1)1' + 2 |G(t)| mt, 1)1 |at0,(i, 1)1 ] 

(2.14) +2|C(t)| |19t0,(t,-)||-i,,. |P,(/3(t,-)e„(t,.))|i 

+ 2^(t)Im(«'4t,-),at0,,(i,-)) Vie [o,r], 

where 5* :— inf[Q S > 0. In order to bound properly the quantity \\dtdh\\-i,h, first use (|2.7p to obtain 

(2.15) ll**(i,-)ll< Ch''+^[\\u{t,-)\\r+i + \\dtu{t,-)\\r+i] yte[0,T]. 
Then, use of ([22]) and (^^5)) in gives 

|(at^/.(t,-),'/')| < |a(t)| [s{t) \^l{t)\ \dteh{t,i)\ + {\G{t)\ \^,{t)\ + 1) |0,(t,.)|i] |0|i 
+ |/3(^,•)U|j^^,(^,•)||||0|| 

+ Gr+i(|i9tw(t,-)||r.+i + h(t,-)llr+i) II0II V0e5,„ V<e[0,r], 



which, along with (|2.ip , yields that 

2 1^(01 \\dt9h{t, < 2 S{t) \^i{t)\ \dte^{t, 1)1 

^^■^^^ + C [ \eh{t, oil + h' + '{\\dtuit, ■)\\r+l + \\u{t, ■)\\r+l)] Vie [0, T]. 

Thus, combining (|2?T4l) . ([2?T6l) . jSH]), jMl), and ([2^ . we arrive at 

^|e/.|i< C [|0,|? + /i2('^+i)(||a,u||2+i + ||?.||^+i)] +2Um(vI..,9t0;,) on [0,r]. 
Since 0/i(O, •) — 0, integrating with respect to t in the inequality above yields 



\Oh{t,-)\i<C 



t 

2 





Im 



{2eW(*.(t,-),0/.(t,-))-2 eis){^.is,-),9his,-))ds 
as){dt^As,-),Oh{sr))ds} ViG [0,T]. 



Using in the above the Cauchy-Schwarz inequality, (|2.ip . and (|2.15|) . we obtain 

(2.17) \Oh{t,-)\l<c I \eh{s,-)\\ds + ch^^-+^^Vj^{t) v<G[o,r]. 



The estimate (|2.12l) follows from (I2.17P using Gronwall's lemma and (|2.ip . □ 

A simple consequence of this superconvergence estimate and the approximation property (|2.7p of the 
elliptic projection is the following convergence result: 

Theorem 2.5. Let u be the solution of (Af) and Uh its semidiscrete approximation defined by (|2.10p - (|2.1ip . 
Assume that fi(t) < and S{t) > for t G [0,T]. Then 

( .M|2 



(2.18) \\uf,it,-)~uit,-)\\+h\\uhit,-)-uit,-)\\i< Ch''+' ^\\u{t,-)\\Ui + ]^ ^M{T)dT^ ViG[0,r], 

where is the function defined in the statement of Proposition \2.4\ O 

Therefore, taking into account the relation of a, /i and S to the function s(t) describing the bottom 
topography, we conclude that the error estimate of Theorem 12.51 holds in the case of domains with upsloping 
bottom profiles, i.e., when s{t) < for t G [0,T]. 

Remark 2.1. The superconvergence estimate (j2.12p . (j2.2p . and a standard L°° estimate for the error 
of the elliptic projection ([25j) yield as usual an optimal-order estimate of the error \u — Uh\<x> on [0,T] (c/. 

ESI). 

2.2.2. Crank- Nicols on fully discrete approximations. Let iV G N and (i")^^o t>e the nodes of the partition 
of [0, T] where, = 0, = T and < for ti = 0, . . . , iV - 1. Define fc„ - for n = 1, . . . , iV, 

^„+i _ t"+t"+' for n = 0, . . . , iV - 1, and fc maxi<„<„ fc„. We set u" u(t", •) for n = 0, . . . , iV, 
where u is the solution of [N). Finally, for sequences {V"^)^^q, we define dV"^ :— -^{V"^ — \/™^i) and 
ylT/™ = 1 (T/™ + y"-i) for TO = 1, . . . , M. 

For n = 0, . . . ,N, the Crank-Nicolson method yields an approximation UJ^ G Sh of M(t", ■) as follows: 

Step 1: Set 

(2.19) := Rhuo. 
Step 2: For n = 1, . . . , N, find C/,^ G Sh such that 

(9[/;;\x)=ia"-^A^"-5 [5"-^5C/^(l) + G"-5^[/;;^(l)]^ 

-ia"-^B{AUj:,x)+i{l3''-^AUj:,x) + {f''-Kx) Vx e 5,, 
where S'"-^ ;= ^(t"-^), ^"-^ :^ ;u(t"-5), a"-^ := a(<"-5), G"-3 := G(i"-5), := /(<"-5,.) and 

:=/3(t"-3,.). 



We first examine the problem of existence and uniqueness of the fully discrete approximation 



Proposition 2.6. Let n e {!,..., A^} and suppose that UJ^ ^ ^ Sh is well defined. If S^^ 2 > and 

< 0, then, there exists a constant C„ such that if < Cn, then UJ^ is well defined by (j2.20p . 

Proof. Since (|2.20p is equivalent to a linear system of algebraic equations with unknowns the coefficients 
of with respect to a basis of Sh, existence and uniqueness of C/^ will follow if we show that if there is a 
V G Sh such that 



(2.21) 



then V ~ 0. Set — -^V in (I2.2ip . and then take imaginary parts and use the arithmetic-geometric mean 
inequality and (|2.9p to obtain 

1 Re(G"-5)|i/(l)|2 



\V\i =Ai""2 fc„ 



(2.22) 



2 5"~2 
-S"-i 



y(i) 



fen 



<|M"-'|fc: 

For G 5/1, we use ([23T|) . ((2?2)) and ((27T|) to obtain 

|(£,0)| <|a"-^||M"~^l [5"-^ 

which yields 

(2.23) II 
where ^ [ |a"-^ | + C |/3' 

(2.24) |l/|2(l-A:„ 



jGri^|y(i)|2 



|a"-7| 



+ iG"-^|iini)i] 1011 
i^iii-^ii, 



£||-V^< K-^|K-^|5"-5 



V^l) 
fe„ 



— I ' — r- 



G.lV^li, 

5 1 |^"-^| |G"-^|]. Using arm and (l^?^ . (1^ gives 

__ 1 

|i,oo I I I c"—^ I /^n-^ |2 



VI/^"~^l^"-^l/3"-^lt, 



<0, 



which ends the proof. □ 

In particular, if we suppose that (3 is in G([0,T], W^'°°{D)), that a, /i, 5, G are continuous functions 
on [0,T], and that S{t) > and fj,{t) < for t e [0,T], (i.e. the upsloping case), then the existence and 
uniqueness of the fully discrete approximation UJ^ follows if fc„ < G, where G is a constant independent of n. 
This follows from Proposition [511] and the fact that the quantity multiplying fc„ in (|2.24p may be uniformly 
bounded with respect to n. 

In the case of a general bottom topography we have: 

Proposition 2.7. Let n e {!,..., iV} and suppose that U^~^ £ Sh is well defined. Then, there exist 
constants Cn,i and Cn,2 such that if ^ < Cn,i and fc„ < G„^2, then UJ^ is well defined by (I2.20p . 

Proof. Let dim5;i = J and be a basis of Sh consisting of real- valued functions. It is easily 

seen that existence and uniqueness of UJ^ is equivalent to the invertibility of a matrix M S C"''^ '' defined 
by Mij := M{(t)j,(j)t) for j, i = 1,...,J, where M : Sh x Sh ^ C i s giv en by M{x,(t>) (x,^) - 
ia"-^/^"-5 5"^x(l)0(T) + ^ [ - i(/3"~^X>)-iA'""^a""^G"-^x(l)^ + ia""^S(x>)] for 
Sh- If X e Ker(M) we have Re[M{(f>^, (/)*)] = with 0. Y.j=i 't^i- Then, using (jO)) and (P3)) . we get 

l/3"~^|oo II0.IP + 2|m"-^| |a"-^| |G"-^| 110.11 I0.I1 



< 



l/3"-^|oo +f Im"-^IK-'I|G"-^| 110.11= 



which, under our hypotheses, yields a; = and ends the proof. □ 
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Hence, if /3 G C{[0,T], L°°{D)) and a, /i, G are continuous on [0,T] (i.e. in the case of general bottom 
topography), the existence and uniqueness of Uh follows if we take fc„ < Ci and < C2, for some constants 
Ci and C2 independent of n. 

We next establish the consistency of our fully discrete scheme in the t variable. 

Proposition 2.8. Let u he the solution of (Af). For n — 1, . . . , N , define a" : D ^ C by 

(2.25) =ia"~5U:,:,(t""5,-) + i/3"~^^M" +/"~^ 

Then, 

(2.26) \\a^< C{Kf B'^{u), n=l,...,7V, 

(2.27) -a"\\<C[ (fc„)' + \kn+i - fc„| ] (fc„ + fc„+i) B^u), n = 1, . . . , TV - 1, 

where -B"(u) := ^^^2 ni^X[t"-i.t"] a^irf ^2 (") X]l=2 ™^^[t""i-t"+il ll^^""!!- 

Proof. It follows easily by using the partial differential equation and Taylor's formula. □ 

We prove now that the following superconvergence estimate holds in the fully discrete case. 

Proposition 2.9. Let u be the solution of (M) and (C/^)^=o be the fully discrete approximations that the 
method ()2.19p - (|2.20p produces. Assume that fi{t) < and S{t) > for t e [0,r]. In addition, assume that 
there exists a constant C > such that 



(2.28) 



\kn+i - kn\ < C max {fc^, kl_^_^} , n = 1, . . . ,N - 1 



Then, there exists a constant C\ such that, if max (fc„ Ci) < h, there exists a constant C > such that 



(2.29) 



l<n<N 

max \\UJ: - RhU^'Wi < C {k^ + h''+^) S^(u) V/iG (0,/i.], 

l<n<N 



where ^^{u) I]^=o ™^^[o,ti + max[(,.r] H^i'ulli + maX[o,r] \\dfu\\ + maxte[o,r] \dfu{t, 1)|. 

Proof Let h e (0, h^, := - for n = 0, . . . , iV, C := ^ and ^""^ for n = 1, . . . , iV. 

We use (11201), (1^^ . (Ell), and dM]), to obtain 



(a0;j,x)=ia"-5M 



(2.30) 



where 



5"-^ dei{\) + G"-5 ^0;:(i) - £3" J x(i) 

- ia"-t B{ABl,x) + i (/^h(/3"-^ ^^^;:), x) 
+ (fr-^",x) + i«"-^S(£2",x) Vxe^h, n = l,. 



El := au" _ i?;,(aM") - iP,,[/3"-5(ylu" - i?^(ylw"))], 
£2" :=w(r-5)_^y", 

£3" :=5"-^ [atu(<"-^l) -<9u"(l)] +G"-5 [it(t"-5,i)_^w"(i)]. 
Using Taylor's formula and (|2.7p . we deduce the following estimates: 



(2.31) 

(2.32) 

and 

(2.33) 



< Ch''+^ 



dtu{s, ■) ds 



max u ,.+1 



+ C|/3"-5|oo 

r+l 

max ||9tu||r+i 



\Au^\ 



r+1 



\£^\<Ckl 



< Ckl max 



max 192^(^,1)1+ max \d^Mt^)\ 

te[t"-i,t"] te[t"-i,t"] 
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for n = 1, . . . , N. Set x — dO^ in (|2.30p . and then take imaginary parts to obtain 

(2.34) + fc„ m{l)\' + 2 |G"-^| \A6l{l)\ \dei{l)\ + 2 jfjl ] 



,...,7V, 



where S"* := inf^p^T] -S"- 

Now let us estimate \dd'^\\-x.h. For e 5"^, (|230| - ((233| . ([2:26)) . Jl^), and (El]) give 

1(90^ ^)| < K-^l 5"-^ Ml 

+ c|^0;:ii|^|i + c(;i'-+i + fc2) 1^1^ n = i,...,7v, 

where := maX[o,r] ||u||r+i + maX[o,T] ||9fiij|r+i + maX[o,r] i + maxjo.Tj + maxtejo.Tj 1)|. 

Hence, we conclude that 



(2.35) 



2fc„|r~3|||a0;j||_i,,. <2A:„|Ai"-^|5"-2 \ddl{\)\ 

^CK\AQl\x + Ckn{lf^^ ^kl) Si(u), n= l,...,iV. 
Now, combining (|2.35|) and (|2.34p we have 

< K''\i + ckn \Aei\\ + cfc4 (fc„)^ + {h^+' + kl) \Aei\r ] Si(u) 

+ 2 /c„ Re[S(£2", 96i;j)] + 2 A:„ l^^i - ct", 96i;j), n = 1, . . . , TV, 
from which there follows that for some constant Ci > 

(1 - Ci fc„) < (1 + Ci fc„) + C2 A:„ (r+i + fc2)2 (Si(^i))2 

(2.36) 

+ 2 /c„ Re[S(£:2", ae;^)] + 2 A:„ ^ im(f - a", 96*;^), ?i 1, . . . , TV. 
To continue, we assume that max (Ci fc„) < ^ which allows us to conclude that < g'^'-^i'^" for 

l<n<iv^ ' - 3' 1-Cifc„ - 

n = l,...,N. Hence, (|236)) yields 

qn|2 ^ 3Cik„ \nn-l,2 . C2 k„ fur+l 1 l,2\2 /„,^^2 



2k„ 



1-Ci fe„ 



n = 1, 



,N. 



Next, we define A^* := °''p(^^'^^^g^J+l ^0 g^j-^j simple induction argument to arrive at 



K\l<c2iE^{u)rY.k,\-{h^+' + k^y 



which yields 

(2.37) 

where 



2 ^ k, X] [Re[B{£^,dei)] + Im(£f - a^ 



m< C{h^+' + k'riE,iu)y+T: + Tl\ n=l 



T;:=2^A;;Re[S(£^",0i-0r')], 

n 
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n = 1 



,...,N, 



, . . . , J. T , 



First we observe that 



(2.38) U 



n-l 



for n = 1, . . . ,N. Since 

- < C {k, + k,+i) [ {kjf + |fc,+i - fc,| ] J = 1, . . . , TV - 1, 

with S2(w) := niaX[o,T] \dtu\i + maX[o,T] we see that (j2.38p . (|2.32p . and (|2.28p yield 

(2.39) |T;| < 22(7.) max l^rii, n^l,...,N. 

l<m<n 

In addition, we have 

n-l 



(2.40) 



n—l / n \ 

25: exp 3Ci ^ fc. [ exp(3C.^^,^.)-i+C./c., lm(£r'-a^+\0i) 

j=l V ^=J+2 / 



n-l 



+ 2 5^(e-^ - e+^) AJ+i Im(f - 0i), n = 1, . . . , iV. 
Observing that 

- £i^'\i < C (k, + fc,+i) h-+' Esiu), J = 1, . . . , iV - 1, 

with Es{u) := maX[o,Tj ||5tu||r+i + maX[o,Tj ||(9t^u||r+i, we see that (IOT|) . (PT^ - iP?^ and, (EH]) yield 

(2.41) |r^| < C(fc2 + r+i)S4(u) max n=l,...,7V, 

l<m<n 



where ^4(u) := X]^=o ^^-'^P-ri + ^£=3 ni^^io.r] Now, from (|2.37p . (|2.39p . and (|2.4ip there 

follows that 

K\l< C{h^+' + kY{E,{u)f + C{k' + h^+'){E2{u) + E,{u)) max l^rii, n=l,...,iV, 

l<m<n 

which easily yields 

(2.42) max \e]l\l < C ih-+' + k^ f (E^iu) + E2iu) + E^iujf . 

0,<n<N 

The desired estimate (|2.29l) is then a simple consequence of (|2.42p and (|2.ip . □ 
Now we are ready to prove error estimates in the and norms. 

Theorem 2.10. Let u he the solution of (A/") and {UJ^)!^^q be the fully discrete approximations that the 
method ^J^-^^ produces. Assume that fi{t) < 0, S{t) > 0, for t £ [0,T], that (P?^ holds and 
maxi<„<„(Ci kn) < ^, where Ci is the constant specified in Provosition \2. 9[ Then 

max ||J7^"-7i"||f < C(fc2 + /i'-+i-n S^(u), yhe{0,h^], 

0<n<N 



for ^ = 0, 1, where :Zj^(u) was specified in Proposition \2.9\ 

Proof. It is a simple consequence of (j2.29p and (|2.7p . □ 

12 



We conclude that in the case of upsloping bottoms, the fuUy discrete Crank-Nicolson-Galerkin method 
(|2.19p - (|2.20p yields fully discrete approximations U]^ that converge to the solution u of (Af) at optimal rates 
in the and norms. 



2.3. The Abrahamsson-Kreiss boundary condition. We consider now the (PE) with the Abrahamsson- 
Kreiss bottom boundary condition, i.e. the ibvp (|1.6p . (|1.8|) . (II. 9p . (I1.12p . which we rewrite here, in slightly 
more general form, for the convenience of the reader. For T > given, seek a function u : [Q,T] x D ^ C 
satisfying 

ut = i a{t) u^^ + i (3{t, x)u + fit, x) V(t, x) e [0, T] x D, 
u{t,0) = Vi G [0,T], 
^^'^^ uAtA) = yte [0,T], 

u{0,x) = uq{x) Va; G D. 

We assume again that a : [0,T] IR\{0}, /3, / : [0,T] x D ^ C, uq : D ^ C are given functions. We shall 
assume that the solution of (AlC) exists uniquely and that the data and the solution of (AK.) are smooth 
enough for the purposes of the error estimation. We note that (AlC) may be considered as a special case of 
(Af) obtained by setting fj. equal to zero in (M) . (This does not imply of course that we assume that s is 
zero. We recall that in the Abrahamsson-Kreiss formulation the effect of variable bottom enters explicitly in 
the definition of a and /3, cf. (|1.7p . and in the change-of- variable formula p.Sp .) All the error estimates for 
(AJC) that follow may then be considered as special cases of the analogous estimates in the two preceding 
paragraphs but with some important simplifications. For the convenience of the reader we shall restate the 
results but not prove them in detail; we shall just point out some differences between them and the analogous 
estimates for the problem (Af). It will be seen that the finite element approximations of (AlC) exist and 
satisfy optimal-order error estimates under no further assumptions (except smoothness) on the shape of the 
bottom. 

2.3.1. Semidiscrete approximation. Using the finite element subspace and the notation established in 
paragraph [^m we define the semidiscrete approximation Ufi of the solution of (AlC) as the map Uh ■ [0,T] ^ 
Sh satisfying 

(2.43) {dtUh{t,-),q})^-ia{t)Biuh{t,-),^)+i{P{t,-)uh{t,-),q}) + if{t,-),^) V(^g5,„ ViG[0,T], 
and 

(2.44) UhiO,-)^ul 

where G Sh is an approximation of uq, which may be taken, for example, as PhUq or Rhuq. 

Proposition 2.11. The problem (|2.43p - (l2.44p admits a unique solution in C^{[0,T], Sh)- // / = and 
l3j = 0, then the solution preserves the L^{D) norm, i.e., 

(2.45) \\uh{t,-)\\ = \\ul\\ VtG[0,r]. 

Proof. The first part follows from Proposition 12.31 for /i = 0. The conservation of the L^{D) norm follows 
by taking cj) — Uh m (|2.43p and then real parts. □ 

Theorem 2.12. Let u he the solution of {AK,) and Uh its semidiscrete approximation defined by (|2.43p - 
((2:441) . Then 



(2.46) \\u{t,-)-Uhit,-)\\i < C 



VtG [o,r], 



for e = 0,1 and h G (0, /i^], where r^K;.o(T) := Y.ln=o W^ruir, Ollr+i ^Aic.ii'r) ^jsrir), where 
the function defined in the statement of Proposition \2.4\ 
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Proof. Let h G (0, h^,]. Defining as usual 9h := Uh — RhU we obtain 

(2.47) {dtehit,-),q>) = -ia{t)B{0hit,-),(l))+iiKt,-)Ohit,-),^) + {'i'*{t,-),4') V(^e5,„ Vte[0,T], 

where was defined in tlie course of the proof of Proposition 12.41 Taking (p ~ 9fi in (|2.47|) and then real 
parts, we may prove (|2.46p with £ = in a straightforward manner. The proof of (I2.46P with £ ~ I follows 
the steps of the proof of Proposition if we take = dtOh in (|2.47p and then imaginary parts. □ 

Remark 2.2. Hence, if is taken equal to Phu'^ or RhU^ , Theorem \2.1S\ yields optimal- order estimates of 
the error u — Uh in the or norm, respectively. Also, we note that to obtain the estimate (|2.46p with 
£ = 0, we do not need the inverse inequality (j2.5p . 

2.3.2. Crank- Nicols on fully discrete approximations. We now proceed to the full discretization of (A)C) by 
discretizing the initial- value problem (|2.43p - (|2.44p in t using the Crank-Nicolson scheme. With notation 
introduced in paragraph 12.2.21 we define for n = 0, . . . , approximations Uf^ £ Sh of u(t", ■), the solution 
of {AJC), as follows: 
Step 1: Set 

(2.48) t/o ul 
Step 2: Porn = 1, . . . , iV, find C/^ G 5"^ such that 

(2.49) (9C/^\x) = -ia"-^6(^(7,",x)+i(/3"-^^C/,?,x) + (.r-',x) VxG^h- 

Proposition 2.13. Let n e {1, . . . , iV} and suppose UJ^^^ is well defined. Then, there exists a constant C 
independent of n such that if fc„ < C , UJl is well defined by (|2.49p . Moreover, if f = and f3i = 0, then 

(2.50) \\Un = \Kl n = 0,...,N. 

Proof. Since (|2.49p is equivalent to a dim(S'/i) x dim(S'/i) linear system of algebraic equations, existence and 
uniqueness of UJ^ will follow if we show that if there is a. V £ Sh such that 

(2.51) ^(T/,0) = -ia"-^6(y,0) + i(/3"-^l/,0) Vc^eSh, 

then V — 0. This fact follows easily for fc„ sufficiently small, if we put (f) — V in (|2.5ip and take real parts. 
The conservation property (|2.50p follows from (|2.49p if we select x = -^U'h and take real parts. □ 

Theorem 2.14. Let u be the solution of [AJC) and {UJ^)n=o be the fully discrete approximations produced 
by ^Ml-^M- Then, if ma,xi<|j<jv is sufficiently sttkiII, we have 

(2.52) max \\UJi ~ < C [\\ul - RhUo\\ + (e + h^+^)E^^,„{u)] V/ie(0,M, 

0<n<JV 

wfcere S^k;,o('w) := J2L=o^^^io.t] ll<9r"llr-+i + Tfm=2^^^^[o.T] Wd^uW + maX[o,Tj \\d'^dlu\\. Also, if (HIS]) and 
(I2.28P hold, and maxi<„<„ fc„ is sufficiently small, then 

(2.53) max ||(7;^ - < C [ |K - i?,.uo||i + (fc' + /i") S^;c.i(w) ] V/^e(0,/i.], 

0<n<iv 

where ^Aic.a{u) + maX[o.T] \\dtu\\r+i + maX[o,T] \\dfu\\ + maxjo.T] \\dfdlu\\. 

Proof. First, we modify the consistency argument of Proposition [2T8] defining, ior n — 1, . . . , N , a" : D ^ C 
by = ia"-^^(M:r^(t",-))+i/3""^^u" + /""5+cr". Then, we set 6i;^' UJ^-Rhu" for n = 0, . . . , iV, 

to obtain (j2.30p simplified by setting /i"^^ — and £2 = 0- To obtain (|2.52p we put X = ^-nd then take 
real parts. To obtain (|2.53p we proceed along the lines of the proof of Proposition [2?9l appropriatelv simplified. 
□ 
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3. Numerical experiments 



In this section we present the results of some numerical experiments that we performed using the fully 
discrete Galerkin-finite element methods, defined and analyzed in the previous section, to solve the ibvp 
for the (PE) in domains of variable bottom topography with Neumann and Abrahamsson-Kreiss boundary 
conditions. We also make, in paragraph 3.3, a theoretical excursion with the aim of explaining some experi- 
mental observations made in paragraph 3.2. Recall that in the case of the Neumann boundary condition, i.e. 
for the problem (Af), our convergence results were rigorously established in the case of upsloping bottoms, 
that is when s{t) < for all t g [0, T]. One of our goals in this section is to study numerically the behavior 
of the Neumann boundary condition in the presence of downsloping bottoms and compare the solution of 
(A^) with that of (AJC), for which rigorous convergence results hold for any smooth s{t). In the numerical 
experiments the finite element subspace Sh consisted of continuous, piecewise linear functions defined on a 
uniform mesh, while the temporal discretization was effected with uniform time step. All computations were 
performed using double precision fortran 77. 

3.1. Order of convergence. To test numerically the order of convergence of the fully discrete Crank- 
Nicolson-Galerkin finite element method (henceforth referred to as (FE)) in the case of the ibvp (A/"), we 
took T — I and considered three cases of bottom profiles, namely: 

Case 1: s{t) — —OM + 0.7 (upsloping). 

Case 2: s(t) = OAt + 0.3 (downsloping). 

Case 3: s{t) = 0.2cos(47rt) + 0.2sin(47rt) + 0.7 (oscillatory). 

In (A/") we took a — l/(2s^), /3{t,x) — xt + i{3x + t^), uq{x) = —x{x — 1)'^. The bottom boundary 

condition had the form u^{t, 1) = ^i{t) [Sit)ut{t, 1) + Git)u{t, 1)] + fi{t), where fj.{t) = S{t) = i/Jlly2 , 

G{t) — i{S{t)6{t) — s'^{t)), 6 = ss/2. The nonhomogeneous terms / and /i were chosen so that the exact 
solution of the problem was given by u{t, x) — ~x(x — lY +sin(t)x. To compare the exact with the numerical 
solution we calculated the I2 error at the nodes Xj at T = 1 (taking k = h). Table 4.1 shows the rates of 
convergence of the numerical solution in the three cases. The rate is clearly two in the upsloping case 
(as predicted by the theory), approaches two in the downsloping and seems not to have stabilized in the 
oscillatory case. On the other hand, as predicted by the convergence theory, (FE) when applied to {AJC) 
gave clear second-order convergence. 



h 


Case 1 


Case 2 


Case 3 


1/100 


1.998 


1.638 


1.766 


1/200 


1.999 


1.659 


1.085 


1/400 


1.999 


2.001 


1.556 


1/800 


2.000 


2.012 


2.615 



Table 4.1. Orders of convergence of (FE) for [M] in I2 in three cases of bottom topography. 

3.2. Comparison of {M) and {AJC): The upsloping and downsloping wedge. We first consider the 
ASA upsloping wedge underwater acoustic test problem, see [15| , with rigid bottom given in the original 
variables r, z by the function l{r) = 200 — 0.05r m for < r < 3339 m. The source, of frequency /o = 25 Hz, 

was placed at Zs — 100m and modelled by the initial value iPq{z) — y^^{exp(— (z — Zs)^-^) — exp(— (z + 

2s) x)}' < z < l{0). The water was assumed to have constant sound speed equal to c = Cq = 1500 m/sec 
and no attenuation. In (PN) gB{r) was taken equal to ifco- The problem was transformed by the change 
of variables (|1.6p to an equivalent one on the horizontal strip 0<x< l,0<i<T, and it was solved 
numerically by (FE) in both the {J\f) and {AK.) formulations with h = 1/1000, k = T/1000, T = 3339. (In 
the figures that follow we present the numerical results after transforming them back to the original r, z 
variables. Specifically, we present graphs of the numerically computed field ip, represented as is customary 
in underwater acoustics, by the transmission loss function TL = —20\ogiQ{\^p{z, r)\) + lOlog^Q r dB depicted 
as a function of r at certain depths z.) For this upsloping example we show in Figure 2 the transmission 
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loss curves as functions of r e [0, 2200 m] at a depth oi z = 90 m for both the (Af) and (AlC) models, which 
evidently agree very well. 

depth-90m 


10 
20 
30 
40 

^ 50 

60 
70 
80 
90 
100 

200 400 600 800 1000 1200 1400 1600 1800 2000 2200 

range in meters 

Figure 2. Upsloping ASA wedge; TL as a function of r at depth z — 90 m, comparison of (TV) and (AK.). 

We then considered the analogous downsloping wedge given by l{r) — 33.05 + O.OSr for < r < 3339 m. 
The source, of frequency 25 Hz, was placed at Zs = 25m and modelled as in the upsloping case. In this case, 
we found that the (FE) numerical solution of the problem (Af) apparently exhibited numerical instabilities 
and did not seem to converge as the discretization parameters became smaller. For example, in Figure 3 
we superimpose the TL curves at depth z = 25m corresponding to the (Af) model solved by (FE) with 
h = 1/100, k = T/lOO and h = 1/1000, k = T/1000, T = 3339, with the analogous results obtained by 
(AlC) solved by (FE) with smaller h and k. The (AlC) model, when discretized by (FE), yields reasonable 
results that converge to the solution shown in Figure 4 with dotted line. To make sure that the numerical 
method used for (Af) was not the culprit, we repeated the numerical experiment using a Crank-Nicolson /inzte 
difference discretization for (Af), and found results identical to those of the (FE). We tentatively conclude, 
therefore, that in this realistic downsloping bottom case, the model (Af) allows the growth of instabilities, 
in agreement with the remarks of Abrahamsson and Kreiss in [T] and [2, . 

depth— 25m 

umann h=l/100,k=T/100 I 
umann h=l/1000,k=T/1000 
K h=l/1000,k=T/1000 



500 1000 1500 2000 2500 3000 3500 

range in meters 

Figure 3. Downsloping ASA wedge; TL as a function of r at depth z — 25m. (FE) solutions for the (Af) and (AK.) 

models. 

To check the validity of the (AK.) solution of this problem we compared the results of Figure 3 with those of 
yet another numerical method, the Crank-Nicolson type finite difference code IFD for the (FE), [16], |17j . 
|18j . which has been widely used in underwater acoustic numerical simulations. 
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(b) 

depth=25m 




100 I ' ' ' ' ' ' 1 

500 1000 1500 2000 2500 3000 3500 

range (km) 

(c) 

Figure 4. Downsloping ASA wedge; TL as a function of r at a depth z = 25m. Comparison of (A/") and (AlC), 
discretized by (FE), and IFD with rigid bottom b.c. (a): fo = 25 Hz, (b): Magnification of (a) for r € [1210, 1430], 

(c): /,, = 80 Hz. 

We chose the option of the rigid bottom boundary condition in IFD and solved the problem using Az = 
3.31 m, Ar = 0.17 m, values by which the IFD solution had converged. (The IFD code solves the problem in 

the original r, z wcdgc-shapcd domain). Figure 4(a) shows the superimposed TL curves obtained at z = 25m 
by the (A/") and (AK.) models solved by (FE) with h = 1/1000, k = T/1000, T = 3339 (as in Figure 3) 
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and for the IFD with the rigid bottom boundary condition. The results of {AlC) and IFD agree well. In 
fact, they differ by about half a dB as inspection of a typical window of Figure 4(a), shown in Figure 4(b), 
reveals. (It is worthwhile to note that at a higher frequency fo — 80 Hz the results of (FE)-(A/') approach 
those of (FE)-(yl/C) and IFD, see Figure 4(c)). 

To explain this result we looked closely at how IFD implements the rigid bottom boundary condition and 
found that it does not actually discretize (PN); instead, it uses a different boundary condition obtained by 
replacing the tpr term in (PN) by 2^ipzz + ^{r]'^ — l)^' using the (PE), and then discretizing the ^pzz term 
at the bottom with one-sided finite differences from the interior of the domain. In the next paragraph we 
offer an explanation why this rigid bottom boundary condition yields a stable problem for any monotone 
bottom profile. 

Our tentative conclusion, then, from this experiment is that in the case of realistic, downsloping en- 
vironments, (AJC) and the rigid bottom boundary condition model implemented by IFD apparently yield 
correct results, while the Neumann bottom boundary condition used in (A/"), which retains the term ipr at 
the bottom, allows the growth of instabilities. 

3.3. Using the p.d.e. in the dynamical boundary condition. Let w = w{t,y) be defined for < y < 
s{t), < t < T, and satisfy (jl.H) - (|1.4|) . Replace the term Wt{t, s{t)) in p.4|) by its value given by the p.d.e. 
in p.ip to obtain 

(3.1) Wy{t, s{t)) - s{t) { i wyy{t, s (t)) + [i {t , s (t)) + g (t)] w{t, s{t)) } = Vt G [0, T]. 

In the IFD code, the rigid bottom boundary condition used is a finite difference discretization of (|3.ip . 

To avoid the presence of the second derivative Wyy{t, s{t)) in the boundary condition (|3.ip . we differentiate 
p.ip with respect to y and put p(t, y) — Wy{t, y). (Note that w{t, y) — Jq pit, ^) since it;(i, 0) = 0.) Then, 
the ibvp ITTIl-dOl), im]) becomes 

pt^ ■^Pyy + n{t,y)p + ny{t,y)w yye[o,s{t)], \/te[o,T], 
Py{t,Q) = Q vte[o,r], 

p(t,s(i))-s(0{5Pi/(^'*W) + [i7(i,5(t))+g(i)] w{t,s{t))]^Q Vi G [0,T], 
p{Q,y)^p^{y):=w'^{y) Vye[0,s(0)]. 
(Note that using the (|l.ip at y = and the surface boundary condition w{t, 0) = 0, we obtain that Py{t, 0) ~ 

Wyy{t,Q)^0.) 

In what follows, we shall obtain an a priori bound for the solution of p.2p and then propose a finite 
element method for solving it. With this aim in mind, we perform as usual the range-dependent change 
of depth variable x :— that maps the domain of the problem onto the horizontal strip {{t,x) : t G 

[0,T], X G D}, where D = (0, 1). Consider the transformation 

(3.3) p{t,y) = ^ exp{-C{t,x)) (^p{t, x) ~ (Jt, x) ^ p{t,Od^ 

where the function ( will be specified below. Note that the function 6, defined by 9{t, x) := Jq p{t, £,) for 
{t,x) G [0,T] X D, satisfies the first-order o.d.e. 

9At,x)~C,{t,x)0{t,x) = s{t) eMC{t,x))p{t,xs{t)). 

Solving this differential equation with initial condition 6{t, 0) = yields 



e{t, x) = s{t) exp(C(t, x)) / p{t, C s{t)) de, 
from which we may derive the inverse of the transformation (|3.3p in the form 

p{t,x)^s{t)exp{at,x)) (p{t,xs{t)) + CAt.x) p{t,is{t))d^ 
After some calculations we also obtain that 

(3.4) e{t, x) = exp(C(t, x)) w{t, xs{t)), {t, x) G [0, T] x D. 

18 



Following the ideas of 0], and after analogous computations (see, in particular, (2.7) and (2.8) of t4J), we 
may deduce that p solves a well posed ibvp, in the case of strictly monotone bottoms, i.e. when s{t) is either 
positive or negative for all t e [0,T]. To see this, define first C, as in [4J, by the formula 

(3.5) C,{t,x)^ ^{(j{t)-l)k{t)s{t)x^ V (t, x) e [0, T] X ;D, 

where cr(t) := '^''^t(ty^ ^ + ^(^) > 0, where e is a positive constant, and a[t) := 1, or equivalently C = 0, 
if s(t) < 0. Then, in the transformed domain, and expressed in terms of the new field variables p and 6, the 
ibvp p.2p becomes 



Pt^TiT) + ^(^' ^) P^ + ^) + G{t. x)\p + G,{t, x)e \f {t, x) e [0, T] X D, 

p.{t,0) = V<e[0,T], 
^'■'^ i^p^t, 1) = ^-"-if,^^-'^^ Pit, 1) ^ ^(t, 1) VtMO,T], 

p(0,x) =po(a;) yx<ED, 

where po(a:) - s(0) exp(C(0, x)) K(a:s(0)) + C.(0, x) /; «;[,(^s(0)) d^] , ^(0 = 2s2(t), i?i(t) = i^^ggg^, 
B(t, x) = a: Ig - ^ (t, x) , G{t, x) = Ct x) - a: (t, x) + i 7(t, xs{t)) + 5^ [ (C. (t, x)f^ (^^ {t, x) ] , 

R2{t) = [git) ~ C,tit,l)]Riit) + Cx(i, 1). (Recall that eit,x) = !q pit,Cl d^. In addition, note that (l331l 
yields that £? is real- valued and is given by i?(t, x) = "'(0; that i?(t, 0) = and -Bi:(i, a;) = Bit, 1) = 
^^CT(i). It is easily checked that 1 — i?i (i) i?(t, I) 7^ for t e [0,r].) We may now prove the following result. 

Theorem 3.1. // the bottom is strictly monotone, the ibvp (j3.6p is -stable. 



Proof. Multiply the p.d.e. in p.6p by pit, x), integrate with respect to x in [0, 1], use integration by parts, 
and take real parts to obtain 

+ ReiG,it,-)0it,-),pit,-)) 

+ ^,iB^it,-)pit,-),pit,-)) +ReiGit,-)pit,-),pit,-)) V<e [0,T]. 

Using the Cauchy-Schwarz inequality, the arithmetic-geometric mean inequality, and noting that |0(i, 1)| < 
\\pit, •)||, \\9it, •)|| < \\pit, •)|j, we see from the above that for any ^ > there exists a constant > such 
that 

i\\pit,-)r<{nTit)-^Bit,l)+^) b(t,l)P + Q||p(t,-)|P Vie[0,T]. 

Since -j^^j-pj — ^Bil,t) < for t E [0, T], we may chose ^ sufficiently small to make the first term in the 
right-hand side of the above negative. Hence, by Gronwall's lemma, we conclude that \\pit, •)|| < C||po|| for 
t G [0,T], which ends the proof. □ 

Now, we can define a semiscrete approximation ph : [0, T] Sh of the solution p of problem (|3.6p by 

PhiO,x)^plix) Va;G[0,l], 

and 

idtp,,^) Bip„ 0) + [ '-^iX)'^'' puit, 1) - ^-^n:a^^^'^^ (^r^it, i)] m 

+ ( Bit, ■) d^ph,^) + ( [B^it, •) + Git, ■)]ph,(f>) 
+ iG,it,-)0h,<f>) y(t>^SH, Vte[0,T], 

where 9hit,x) :— Jq Phit,£,) d£, and p° e 5*^, is a given reasonable approximation of po. Consequently, using 
(|3.4p . we see that exp(— (^(t, x))9hit, x) is an approximation of the solution wit, xsit)) of the ibvp (|l.l[) - (|1.4p . 
Also, it follows, as in Theorem 13. 11 that there exists a positive constant C such that \\phit, OH < G \\p^\\ for 

tG [o,r]. 
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3.4. Growth of solutions of (Af) for various bottom shapes. The final set of numerical experiments 
that we report concern the behavior of the size of the solutions of (A/"), as t grows, in the presence of bottom 
profiles of various shapes. Recall that in |T| it was shown that (Af) is well posed if s is strictly monotone, i.e. 
if s{t) > or s{t) < for < t < T. In addition, downsloping bottom profiles were identified for which the 
solution of (Af) grew exponentially with t. (The fact that problems may arise in case s changes sign may be 
expected, in view of the analogous difficulties encounterd in the (real) parabolic case, cf.e.g. [8].) 

The ibvp (Af) was solved numerically with the (FE) method up to T — 1, with /? = / = g = 0, 
uo{x) = —x{x — 1)'^, < a; < 1, with mesh parameters h = k = 1/500, in the case of the eight bottom 
profiles s{t), < t < 1, labeled (a) to (h) and shown in the left-hand icons of the pairs in Figure 5. (In 
all cases depth increases downwards.) The right-hand icon shows the corresponding, numerically computed, 
i2_norm of the solution of (TV) ||u(t, ■)\\ for < t < 1. (Note that ||w(0, ■)\\ = ^ = 0.062994.) The bottom 
profiles are given for < i < 1 by the expressions: (a) s(t) = e*, (b) s{t) = e^*, (c) s(t) = 1 + (i — 0.5)^, 
(d) 3{t) = l-\t- 0.5|3, (e) s{t) = l-{t- 0.5)3, (f) g(^) = 2 - |2i - 1|, (g) s{t) ^ I + [t - 0.5)^ and (h) 
s{t) ^l + t^. 





Figure 5. Behavior of the L^-norm of the numerical solution of (Af) as a function of t for various bottom profiles 

s{t). 
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Only (a) and (b) correspond to strictly monotone profiles for which the theory of [li properly applies. In the 
cases (c), (d), (f) there is a change in monotonicity, in (e) and (g) we have that s{t) ~ s{t) = at t = 1/2, 
while in (h) there holds that s(0) = s{0) = 0. (In the case (f) a i-mesh node was placed a.t t — 0.5, where s 
fails to exist.) 

We observe that the solution maintains a small L^-norm in upsloping, like (b), or eventually upsloping 
bottoms, as in the cases of the trenches (d) and (f). There is a considerable growth of in the examples 
wherein the bottom profile is eventually downsloping, see (a), (c), (g), and (h), in agreement with the 
observations in [T], [5]. We note that in the case (g), an apparent singularity develops at t — 1/2, where 
the bottom curvature changes sign (with horizontal tangent) and the bottom becomes downsloping. This 
apparently causes the L^-norm to grow violently for i > 1/2. A relatively weaker, but sizeable growth is 
also observed in (h) , where the bottom is such that s — s — Oait — and is monotonically downsloping for 
t > 0. One cannot be of course certain about the existence of a singularity at i = 1/2 in the case (g), given 
that the (FE) code does not at present possess an adaptive refinement capability in x and t. However, when 
the experiment was repeated with k = h = 1/800, it was confirmed that the onset of rapid growth occurred 
at about t = 1/2; for this mesh size, became of order 0(10*) at t = 1. 

4. A PARABOLIC PROBLEM WITH A DYNAMICAL BOUNDARY CONDITION 

Here we consider the model one-dimensional (real) parabolic problem ()1.14p with a dynamical boundary 
condition analogous to that of {J\f) , which we re- write here for ease in reading: We seek a real- valued function 
u : [0, T] X [0, 1] M, such that 

ut = a{t) u.,, + f3{t, x)u + fit, x) V [t, x) G [0, T] x [0, 1], 

u(t,o) = o v<e[o,r], 

(4 1) 

^ ^ a{t)u,{t,l) = e{t)ut{t,l) + S{t)u{t,l) + g{t) Vte[0,T], 

u{x,0) — uq{x) Va;G[0, 1], 

where f3 : [0,T] x [0,1] ^ M, / : [0,r] x [0,1] ^ M, 5 : [0, T] ^ R, g : [0,T] ^ R, a : [0,r] ^ (0, -hoc) 
with a* := inf[o,T] a > 0, e : [0,T] ^ R, uq : [0, 1] R, are given smooth functions. We shall construct 
and analyze Galerkin-finite element approximations for the solution of (|4.ip , considering two different cases 
depending on the sign of e. 

4.1. The dissipative case. The dissipative case is characterized by the assumption e{t) < for t G [0,T]; 
the problem is well posed, see e.g. [T^. We assume that its solution is smooth enough for the purposes of the 
error estimates to follow. We adopt the notation and the assumptions of Section [27TJ restricting ourselves to 
the real case, and avoiding the inverse inequality (|2.5p . 

4.1.1. Semidiscrete approximation. Find : [0,T] Sh, a space-discrete approximation of u, requiring 
{dtUhit, •), x) = [e{t) dtUhit, 1) + S{t) Uh{t, 1) + g{t) ] x(l) 

-a{t)B{uh{t,-),x) + {P{tr)uh{t,-),x) + {f{t,-),x) Vx € ^h, Vi€[0,r], 

and 

(4.3) Uh{0,-)^ Rhuoi-). 

Proposition 4.1. If e < 0, then the problem (|4.2|) - (|4.3p admits a unique solution Uh S C^{[0,T]-, Sh)- 

Proof. The result follows if we argue along the lines of the proof of Proposition 12.31 □ 

Proposition 4.2. Let u be the solution of (|4.ip and Uh its semidiscrete approximation defined by (|4.2p - (|4.3p . 

If£{t) <0 forte [0,T], then 

(4.4) \\uh{t,-)~Rhuit,-)\\l< Ch^'^''+^'> (^J^TUT)dT^ Vie[0,T], V/ie(0,/i^], 
where To{t) \\u{t, •)||2_^i + ||9tu(r, Ollr+i- 
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Proof. Let h G (0, /i*] and 0h := Uh - Rnu. Using the p.d.e. in ((^ and ([^ we obtain 

, , {dten{tr),x) = W)dt6u{t,l) + 5{t)eH{t,l)]x{l) - a{t)B{eh{t,-),x) 

(4.5) 

+ W,-)eh(i,-),x) + (*4ir),x) Vxe^ft, Vte[0,T], 
where $^ := [9tM — Rh{dtu)] — (3 {u — Rhu). First we observe that from (|2.7p it foUows that 

(4.6) ll<i>4^,-)ll < [||u(t,.)||r+i + ||atu(t,-)||r+i] Vie[o,r]. 

Setting X = ^'h in (14.511 and using (j2.3p . we obtain 

OIP < ^ |9A(i, 1)P - 2a. |0,(t, •)!? 

+ 2 ||<i>,(t, Oil ||0h(i, Oil + 2 max \(i{t, ■)\ ||0,(t, Of 

[0,1] 

+ 2 (2|<5(0|+e|£(t)|) ||0„(i,OII|e/.(i,Oli Vie[0,T], Ve > 0, 
which, along with (|4.6p . yields 

(4.7) e^||^,(t,.)f < |e(t)||a,0,,(i,l)|2 + C(e + e3)||0,,(t,Oll' + e/i'''^+'^r,(O V^e[0,r], Ve > 0. 
Set X — dtdh in (14. 5p to obtain 

i [a{t) \9,,{t, - 6{t) \9h{t, l)n = -2 |£(i)| \dMt, 1)P - '5(t) |0,.(t, 1)^ 

+ d{t)\9h{t,-)\l^2\\dt0h{t,-)f 

+ 2 (/3(i, Ohit, ■),dt9h{t, ■)) + 2 ($4i, 0, 5t(^,,(i, 0) Vie [0, T], 
which, in view of (|2.3p and (|4.6p . yields that 

(4.8) ^[a(t)|0,(i,OI?-<5(t)|^?h(i,i)n < ~2|e(^)||a,^,(^,l)|2 + c||^?,(^,OII? + /^'('^+''^„(^) vte[o,T]. 

For positive e we define 

(4.9) vS) -.^ e\\eh{t,-W + a{t)\eh{t,-)\\- 5{t)\eh{t,l)\'' V<e [0,T]. 
Then, applying the trace inequality ()2.3p . we have 

ye{t) > e \\e,,it, of + a. \9,{t, 01? - 2|,5(i)| mt, 0|| 

^\e,{t,-)\l+[e-^^) mt,-)r Vte[0,T]. 

If Co := + maX[o,Ti |i5p, (liTTO)) yields that 

(4.11) u,,it)> ^mt,-)\\l vte [0,r]. 

Now, setting e — eo in (|4.7p and then adding the resulting equation with (|4.8p . we obtain 

(4.12) f^iy,Jt)< C,y,Jt) + (e„ + l)h^^^+^^r^{t) yte[0,T]. 

Since 9h{0, ■) = 0, the bound ((44| follows from (|4.12p via Gronwall's lemma and (|4.1ip . □ 
A simple consequence of ()2.7p and (|4.4p is the following optimal-order error estimate. 



(4.10) ^ 



Theorem 4.3. Let u be the solution of ()4.ip and u;i ite semidiscrete approximation defined by (l4.2p - ()4.3 
//e(t) <0 forte [0,T], f/ien 



(4.13) ||u,,(t, - u(t, Oil + ^ \\uh{t, ■) - u{t, Olli <Ch 



r+l 



\Ht,-)\\r+i+(^J^ T^{t) dr 



yte[Q,T], 



where is the function defined in the statement of Proposition \4-'^ □ 



4.1.2. Crank- Nicols on fully discrete approximations. We use the notation of paragraph 2.2.2. For n — 
0, . . . ,N, the Crank- Nicolson method for the problem (|4.ip yields an approximation ?7^' G Sh of u(t", •) as 
follows: 

Step 1: Set 

(4.14) C/O RhUo. 

Step 2: Forn = 1, . . . , iV, find G 5"^ such that 

{dUj:,x)= [e'^-'^dUj:{l)+S--iAUj:il)+g^^-i] x(l) 



(4.15) 



where a"-5 := a(t"-5), 5"-2 := ,5(t"-2), e«-2 e(t"-2) g»-2 .g(i"-5), /"-2 := /(r~2,.) and 
:= /3(t"-5,.). 

Proposition 4.4. Lei n G {1, . . . , A^} anrf suppose that UJl~^ G S*/! is well defined. If e^~^ < 0, then, there 
exists a constant Cn such that if kn < Cn, then UJ^ is well defined by ()4.15p . 

Proof. It is enough to show that if there is a y G S/i such that 



(4.16) 



5"" 2 



V{1) 



(1) 



s(y,0) + i(/3"-^F,</)) v^g5^, 



then V = 0. To arrive at the desired conclusion, first set (j> — V in (|4.16p and use (|2.3p to obtain 

llT^ir + < ^ (2 \V\i - a"-5 |y|2 + H^lp) . 

Then use the arithmetic-geometric mean inequality, to get 

lll^f (1-A:„7„) <0, 

— 

a 2 

The following consistency result is analogous to that of Proposition [ 
Proposition 4.5. Let u he the solution of (|4.ip . For n — 1, . . . , N , define ct" : D ^ M fey 
(4.17) ""X"'' ^a"-^u^^(t"-5,-)+/3"-^yt»" + /"-^ 

T/ien 



1 / 1 I (5"'~ 2 I \ 1 

where = ^ [ l/?""^ |^ -(- ' 1 . This yields y = if we require, for example, kn < jipp- D 



(4.18) 



< C(fc„)2 niax \\dfu\\+ max 



[t"-i,t" 



n = 1, 



,7V. □ 



Proposition 4.6. Let u fee t/ie solution of (|4.ip anc? (?7/")^=o fee t/ie /wZZy discrete approximations that the 
method (|4.14p - (j4.15p produces. Assume that e(t) < /or i G [0, T], and that (I2.28P holds. Then, there exists 
a constant Co > such that: if max (fc„ C^) < ^ there exists a constant C > such that 

l<n<N 



(4.19) 



max \\UZ~Rhu''\\i < C [k^ + T^(u) V/iG 

l<n<N 



where Td(u) := I]l^o ™^X[o,t] -f- X]L2 ^^^[o.r] ll^twlli +Y.'e=2^^^t<^io.T] |9tw(i,l)|- 

Proo/. Let /i G (0, h^, 5" := (5(t"), a" a(t") and 6',^ := L^/^ - for 0,...,N. We use (HTB . 

dHH), (EH), and dSill), to obtain 



(4.20) 



(90^, X) = e"-' 5^(1) + ^"-^ ^^^'(1) - 2J x(l) - S(^(?;:, x) + (/?"-' ^C, x) 



+ (Zr-a",x)+a"-^S(22",x) Vxe5,, 71 = 1, 
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,iV, 



where := ] , Z^' := and Z^' := e"-^ 1)- 

[u(<"-5,l)_y^u"(l)]. Using Taylor's formula and jlJl), we deduce the following estimates: 



(4.21) 
(4.22) 

(4.23) 



< Ch 



\Z3\<Ckl 



r+1 



max ||u||r+i + max ||9tit||r+i 

2 „„„ ia2. 



_te[t" i,f"] 

for n = 1, . . . , A^. The proof now proceeds in four steps 



\Z^h < CK max \dtu\i, 
max _ 1)1 + max \dlu{t,l)\ 



Step I: Set x = -^^h m (I4.20|) and use ()2.3|) and the Cauchy-Schwarz inequality, to obtain 

iicf + ie"i < ii^^r'ip + k"-'i 

+ (e«-^ _ e«) |0;j(l)|2 + (e"-i - e«-3) i(l)|2 
+ 2fc„ [2 W^-'^\\\A6l\\ \A6l\^ - a, \Aei\l + l/?"-^^ M^^jf 
+ V2|Z^'|P0;j||^|^0;j|f 

+ (||a"|| + ||zri|) \\Aei\\ + a'^~'^z-\,\Aei\, 

which, after the use of the arithmetic-geometric mean inequality, yields 



k"IK(i)P< 



(4.24) 



+ \K-'\? 

n\\2 I w ■-rn\\2 



izr 



rn|2 



n = 1, . 



n = 1, 



Step II: Now set x = in to get 



+ 2- 



\f3''-^\^\\A0U\\d9U+a''-^B{Z^^\d0^) 

(l|2ril + lk"ll) \\de]:\\ - zs deui)' 



which, after the use of the arithmetic-geometric mean inequality, yields 

(4.25) +2K \B{z^i,dei)-^z^dei{\ 

a 2 



71 = 1, , 



Step III: For p > 0, we introduce the quantities 

(4.26) V"; := p ( llCir + k"l IC(1)I' ) + K\l - |0r(l)^ = 0, . 

Now, using (|4.26p and (|2.3p we have 

vT>p\\on^ + WH\l~2\^\e^\^\\en 

m — 0, 



>hK\l + \Kf 



n 2 
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1 

Thus, choosing p — pQ :~ ^ -\- 2maX[o,T] ^27 we obtain 

(4.27) VZ>k\K\\l m^O,...,N. 



Step IV: Combining (jiTM)) . g^, ([23]), g^, g^, g^, and glS]) we obtain 



(4.28) 



2 



, n=l,...,7V, 



where Ti(w) := maxjo ||u||r+i + niaX[o ||9fu||r+i +niaX[o.T] ||5j^i*||i+maX[o.T] m|| +maxtg[o \dfu{t,l)\ + 
niaxt£[o,T] \dfu{t,l)\. Using (|4.28p and (|4.27|) we conclude that there exist constants Ci > and C2 > 0, 
such that 



(4.29) 



(1 - Ci fc„) v;; < (1 + Ci k^) v;-i + C2 k (/^'■+' + klf {T^{u)f 



■-\- 2 fc^ 



B{zi^,dei)-^z^dei{i) 



, n=l,...,N. 



To continue, we assume that max (Ci fc„) < i which aUows us to conclude that < e3Cifc„ f^j. 

' — 3' 1-Ci fe„ — 



l<n<w 

71 = 1, . . . , TV. Hence, (|429l) yields 



2 fe„ 



Po 



Letting A^J' 



' l-Ci fc„ 

cxp(3Ci fcjj) 
1-Ci fc,- 



n = l,...,iV. 



a 2 

and using a simple induction argument we arrive at 

n 



— ' ■'I a-' 2 

J = l 



n = 1 



, . . . , J. T , 



which yields 
(4.30) 



V" < C {h^+' + k'Y {Ti{u)Y + Ti" + Ta", 71 = 1 



, . . . , ^ . , 



where T^" 2 Y.''=i ^7 ^i^i'^l ^ and ■= -2 ELi -^a - ^T^Cl))- First, we proceed 

as in bounding the quantity T" in the proof of Proposition 12.91 to get 

(4.31) W'\<Ck^T2{u) max |6l,™|i, n = 1, 

l<m<7i 

where T2(it) := maxjo r] I'^^'^'li + maxjo j,] In addition, we have 

Tl-l 



.,...,7V, 



i=i 



(4.32) 



ra-l 
n-l 



cxp(3Ci 



exp(3d fcj + i) — 1+Ci kj _ Ci fcj + i 
1 — Oi kj 1 — Ci ^j + i 



Observing that 



\Zi - Z^+i| < C {k, + fc,+i) [ {k,f + \k, - fc,+i I ] T3(u), J = 1, . . . , TV ~ 1, 
with T3(m) EL2max,e[o,T] |9|w(i, 1)|, we see that gSaJ, (|2:28l) . g:23l) and ([231) yield 
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(4.33) ITa"! < Cr T3(u) max |16i;nii, n^l,...,N. 

l<7n<n 

Now, from (g^Ql), (|43T|) and (1433]) there follows that 

V;„ < C(r+i + fc2)2(Ti(u))2 + cfc2(T2H + T3(u)) max H^rili, = 

l<m<n 

Use then (|4.27p to arrive at 

max ml<C{h^+'+k'r{T,iu) + T,iu) + Tsiu)f, 

0,<n<N 

which is the desired estimate (|4.19p . □ 

As a simple consequence of (|4.19p and (|2.7p we obtain the following optimal-order error estimates in 
and norms. 

Theorem 4.7. Let u he the solution of (|4.ip and (C/^*)^=o /wZ^?/ discrete approximations that the 

method (|4.14p - (l4.15p produces. Assume that e{t) < for t £ [0, T], that (|2.28p holds and max kji Cjj ^ , 

l<n<N 

where Co is the constant specified in Proposition \4-. 6] Then 

max \\UJi -u"\\i< C{k^ + h''+^-'^) Tn{u), V/ie(0,/i^], 

0<n<JV 

for ^ = 0, 1, where To{u) was specified in Proposition \4-(>\ O 

4.2. The reactive case. In this paragraph, we propose finite element approximations when the dynamical 
boundary condition in (|4.ip is of reactive type, i.e. e{t) > for t e [0, T]. According to [24], [7], the problem 
is well posed only in the one-dimensional case. To construct a finite element method for this problem we 
follow the idea (cf. paragraph 3.3) to replace the term Ut in the dynamical boundary condition using the 



partial differential equation in (|4.ip . Hence we obtain: a{t) Uxx{t, 1) — ^^^jUx{l, t) — + P{t, 1) u{t, 1) — 

+ f{t, 1) for t e [0, T]. Then, to use this as a boundary condition, we formulate a variational formu- 
lation using B{-, •) instead of the L^{D) inner product (•, •). Of course this approach works also if e{t) < 
for t e [0,T]. 

4.2.1. Preliminaries. Let r g N with r > 3, and Sh be a finite-dimensional subspace of M'^{D) consisting of 
functions that are polynomials of degree less or equal to r in each interval of a non-uniform partition of 
D with maximum length h G (0, h^]. It is well-known, [9 , that the following approximation property holds: 

(4.34) ini \\v ~ xh < C h'-^ \\v\Ui, V« e H^+i(Z?), V/i e (0, /i,], s = l,...,r. 

We introduce bihnear forms B*, 7* : H'^{D) x H'^{D) R given by B*{v,w) := {v",w") and j*{v,w) := 
{v",w") + {v',w') for V and w G H'^{D), and set |t>|2 :— for v G H'^{D). Also, we define a new elliptic 
projection i?^ : H^{D) Sh by 

(4.35) j*iRlv,w)=j*iv,x) Vxe5„. 
Lemma 4.8. The elliptic projection i?^ has the following property 

(4.36) iRlvy{l)^v'il) + iRlv^v){l)-^B{R*hV-v,uj) Vvel^iD), 
where lu(x) — x'^ . 

Proof. Let v E M.'^{D) and p = R*^v — v. Since oj E Sh, setting x = in (|4.35p we obtain p"{x) x dx — 
— ^ {p' ,uj'). Then, integrating by parts we get p'{l) = p(l) — p(0) — ^ {p' ,uj'), which is the desired equality, 
since p(0) = 0. □ 
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Proposition 4.9. The elliptic projection has the following approximation properties: 

2 

(4.37) Y.^'\\^hy-^\\i< Ch^+'\\v\U+, 

1=1 

and 

(4.38) \{Rlv ~ vy{l)\ + \{Rlv - ^;)(1)| < C h^ \\v\l+, 
for s = l,...,r, V e H"+i(D) and h E {0,h^]. 

Proof Let h £ (0,/i*], s £ {1, . . . ,r}, v e H^+i(D) and e = R*^v - v. Using (|051) we have 7*(e,e) = 
7'^(e, X — w) for X G which along with (|4.34p yields 

(4.39) |e|2 + |e|i < C/i^-i||«IU+i- 
Now, let w e H^{D) such that 

-w'" + w' ^e' in D, 

(4.40) 

^ ' w{0)^w"{l)^w"{0)^0. 

It is easily seen that (I4.40p conceals a standard two-point boundary-value problem with respect to w' and 
thus existence and uniqueness of its solution follows in a straightforward way; in addition we have that 

(4.41) \\wh < C \e\i. 

Thus, we obtain |le'|p = j*{e, w - x) for x € Sh- Then, we use (gjSl), (jilMjl and ([441]) to get 

|e|? <C (|e|2 + |e|i) h^ \\wh 
<Ch'\\v\\,+i\e\i, 

which yields 

(4.42) \e\i < C h' \\v\\,+i. 

Hence, (I4.37P follows as a simple consequence of (I4.39P and (|4.42p . 
Using gSH), (ESI), (123]), and (ji^T)) we have 

|e'(l)p + |e(l)p<C(|e(l)p-f lle'f ) 
<C|e|? 

which obviously yields (|4.38p . □ 

For later use, we close this section by extending (|2.3p as follows: 

Lemma 4.10. For v G H'^{D) it holds that 

(4.43) \v'{l)f < \v\l + 2\v\,\v\2. 

Proof Let v e H'^{D). Observing that |w'(l)P = J J{v' {xjf x]' dx , we obtain |w'(l)P = Ik'll^ + 
2 J^xv'{x)v"{x) dx, which yields (|4.43p via the Cauchy-Schwarz inequality. □ 

4.2.2. Semidiscrete approximation. We define Uh '■ [0, T] ^ 5^, a space-discrete approximation of u, requiring 

B{dtu,{t, •), X) = { ?g 9,^,(i, !)-[§}+ /3(i, 1) ]u,{t,l)-[^ + fit, 1)] } x'(l) 
(4-44) +f{t,0)x'{0)-a{t)B*{uh{t,-),x) 

+ B{P{t,-)uh{tr),x)+l3{f{tr),x) "^X^S^, Vte[0,r], 

and 

(4.45) M^(0,-) = i?>o(-)- 

Proposition 4.11. If e{t) > /or t G [0,T], t/ien t/ie problem (|4.44p - (|4.45p admits a unique solution 
uh e CH[0,r];5,0. 
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Proof. The result follows if we argue along the lines of the proof of Proposition 12.31 □ 

In the sequel, we assume that the solution of the ibvp (14. ip in the reactive case is sufficiently smooth. 

Theorem 4.12. Let u he the solution of (j4.ip . its semidiscrete approximation defined by (j4.44p - (j4.45p . 

and To be the function specified in Proposition\4-^ V ^i't) > for t G [0, T], then 



(4.46) 



\\uhit,-)~uit,-)\\i < Ch- 



||u(t,.)||.+l+ / To{T)dT 



yt e [o,T], yh e io,h4 



Proof Let h e (0, h^], 9h - R^u, and 77 = R^u - u. Using (|4.44p . (|4.35p and ((i?T|) . we obtain 

BidtOhit, •), x) = { fgy 5=.^'^^, 1) - [I8y + 1)] ^/^(^^ 1) + } x'(l) 
(4-47) - a{t) B*{eh{t, .), x) + B{f3{t, ■) e^{t, ■), x) 

+ B{£^At,-),x)+a{t)B{Rlu{t,-)-u{t,-),x) ^X^Sh, Vte[0,r], 



m 

e(t) 



+ /3(i, 1) r)(t,l). First 



where £rs [^tu - i?^(atu)] - /3 (u - i?;^M) and fH,2(t) := j^dxr]{t,l) 
observe that using (gJll), (|438)) and (gj]), it follows that 

(4.48) |B(£:«,i(i,-),x)+a(OS(ry(t,.),x)l < C/i"^ ( ||w(i, OH.+i + ||9tu(t, OH.+i ) |x|i 
and 

(4.49) \£nAt)x'{l)\ < Ch^\\uit,-)\\r+i\x'{l)\ 

for X e 5/,. and t G [0,T]. Then, set x = ^/i in (|4.47p and use the Cauchy-Schwarz inequahty, (|2?T|) . ()4.48p . 
(lOgi) . SlSl), and (USD, to get 

i ||0/.(i, •)!? < -a. \0h{t, + c [ mt, 01? + r^(0 + |^?h(t, oil mt, OI2] e [0,T], 

which, along the arithmetic-geometric mean inequality, yields 

(4.50) i\9H{tr)\l< C [\eh{t,-)\l + h^^ro{t)] yte[0,T]. 
Since 9fi{Q, ■) = 0, using Gronwall's lemma from (I4.50p we see that 

(4.51) \Oh{t,-)\l< Ch^"- (^J\o{T)dT^ vte[o,r]. 

Finally, we combine (|2.ip . (|4.5ip and (|4.37p to arrive at the error estimate (|4.46p . □ 

4.2.3. Crank- Nicols on fully discrete approximations. For n — 0, . . . , N , the Crank-Nicolson method for the 
problem (|4.ip yields an approximation UJ^ E Sh of u{t'"', ■) as follows: 
Step 1: Set 



(4.52) 

Step 2: For n = 1, . 



, iV, find e 5/> such that 



(4.53) 



B{dUj:,x) 



-^ + r-3(l) 



x'(i) 



+ ^ (0) (0) - a"- 3 (^C/," , X) + S 5 ^C/^ , X) + S (/"- 5 , x) V X e 5^ . 



Proposition 4.13. Let n G {l,...,iV} and suppose that UJ^ ^ G Sh is well defined. If e" 2 > 0, then, 
there exists a constant Cn such that if kn < Cn, then UJ^ is well defined by (I4.53p . 



Proof. It is enough to show that if there is a y e Sh such that 



(4.54) ^B{V,^) 



4^'(1) 
-5 



-4 + ^-^1) 
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y(i) ^'(1) 



for all (j) ^ Sh, then V = Q. To arrive at the desired conclusion, first set (j) = V in (|4.54p and use (|2.fp to 
obtain 



\v\i<'-f 



Then, use (EH), and (ICT) . to get (l-:^c„) < 0, where c„ := |/9"-^|i,oo (f + Cp^) + 2 ( + 



£2 a 2 



-1 \ 2 



/?"- 2 (1) j (7^^ + 1 + iL_|- + ( 1 + iL_|- ) . Thus, assuming that fc„ < 

V = Q.U 



2+c„ 



we easily conclude that 



Theorem AAA. Let u he the solution of (j4.f p and (C/^)^=o f'^^^V discrete approximations that the 

method (I4.52p - (|4.53p produces. If e(t) > for t G [0,r], then, there exists a constant Cfi > such that: if 
max (fc„ Cn) < i there exists a constant C > such that 

l<n<jv 



(4.55) 

where 



max ||?7;;'-m"||i < C(fc2 + r) T«(w), V/ie(0,/iJ, 

l<n<N" 



TJu) := max ||9fu|L+i + max + max \dfd"^u(t,l)\ + ma.x\dfu\2 

[O.Tl [0,T] te[0,T] lO.Tl 



f=0 



m=0 



Proof Let /i G (0, h^, 6^ := UJ^ - and ry" := i?;;^" - u" for n = 0, . . . , TV. Use (l433l) . gTfll and 

dlSni), to obtain 



(4.56) 



I 2 



-/3(i,l) 



(^C)(i) + £;3" + i?4" x'(i) 



- a"- 2 6*(^0^x) + S(r-5 X) 
where cr" : 75 ^ M is defined by l|iT7|) and for n = 1, . . . , TV 

E^' := u{r'^)-Au'\ 

^3 ^-4(5.ry"(l))- f^+/3"-5(l)') ^(^^"(l)), 



,...,iV, 



E: 



■n a 2 



^(u,(fM))-u,(t"-^i) 



*"T + /3"-*(i) ) ^(w(r,i))-u(t"-t,i) 



Using Taylor's formula, (|4.37p and (|4.38p . we derive the following bounds: 

max + max | 



(4.57) B{E'l - a", x) + a"" 5 S(ylry", x) 
(4.58) 
(4.59) 



< c 









Ixb 



|(t"|i < C(A:„)^ max + max 

[t"-i,t"] [t"-i,t"] 



a"-^fi'^(ii;^x) < Cfc2 max l^^b Ixb, 

[f-i.f] 

i?3"x'(l)| < max ||u||,+i|x'(l)|, 



(4.60) 
and 

(4.61) |i?^'x'(l)l < C^A:^ max |a^(i,l)|+ max |a2a^w(t, 1)| ) |x'(l)| 

for n = 1 , . . . , and x ^ ^h- 
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Now, set X — ™ (|4-56P and use (|2.ip . the Cauchy-Schwarz inequality and the estimates (|4.57p . (I4.59|) . 
(|i:BDl) and li^, to obtain 



2^ ( K\l - K-'\l ) < - + C [ \A0'^\l + max I^^IJh 



(4.62) +(fc2 + /i'-)2(T„(u))2 



n = 1, ... ,7V. 



After use of the trace inequahties (|2.3p and (|4.43p . of the inequahty (|2.ip and of the arithmetic-geometric 
mean inequahty, (|4.62p yields the existence of a constant C„ > such that 

(4.63) {l-C^K)K\i< {l + Cnkn)K~^\l + Ckn{k^ + h-f{T^{u))^ n=l,...,N. 

Assuming that max {Cr kn) < \^ and following a discrete Gronwall argument similar to that of Proposi- 

l<n<N 

tion 12.91 we arrive at 



(4.64) max \ei\i < C {k"^ + h'')T a[u). 

0,<n<N 

Thus the desired estimate (|435)) follows easily if we combine (|4?64|) . (|437)) and ((27T|) . □ 
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